Neutron Drops and Skyrme Energy-Density Functionals 



B. S. Pudliner a , A. Smerzi a,fc , J. Carlson c , 
V. R. Pandharipande a , Steven C. Pieper d , D. G. RavenhalP 
a. Department of Physics, University of Illinois, Urbana, IL 61801-3080 
b. Laboratorio Nazionale del Sud, INFN, v Andrea Doria, 95125 Catania, Italy 
Theoretical Division, Los Alamos National Laboratory, Los Alamos New Mexico 87545 
d. Physics Division, Argonne National Laboratory, Argonne, IL 60439-4843 

(February 9, 2008) 



Abstract 

The J 7r =0 + ground state of a drop of 8 neutrons and the lowest 1/2 - and 3/2 _ 
states of 7-neutron drops, all in an external well, are computed accurately 
with variational and Green's function Monte Carlo methods for a Hamilto- 
nian containing the Argonne vi$ two-nucleon and Urbana IX three-nucleon 
potentials. These states are also calculated using Skyrme-type energy-density 
functionals. Commonly used functionals overestimate the central density of 
these drops and the spin-orbit splitting of 7-neutron drops. Improvements in 
the functionals are suggested. 



PACS numbers: 21.10.-k, 21.60.Ka, 97.60.Jd 
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Properties of neutron matter are vitally important in determining the structure of neu- 
tron stars 0, and have a strong bearing on the energies of neutron-rich nuclei, and on the 
r-process in nucleosynthesis 0. It is impossible to extrapolate available data on nuclei to 
the region of neutron matter with sufficient precision using effective interactions. Differ- 
ent effective interactions that fit the energies of laboratory nuclei rather well predict very 
different equations of state for neutron matter Q. In contrast, it appears that modern 
calculations of neutron matter based on realistic models of nuclear forces are much more 
consistent with each other at densities p ^$ 0.16/m~ 3 [|J, and therefore are presumably more 
reliable. The two-nucleon interaction in these realistic models is better determined from the 
scattering data in isospin T=l states than that in T=0, and the uncertainties coming from 
three-nucleon forces and relativistic effects are also much smaller in neutron than in nuclear 
matter. Calculations of uniform neutron matter have provided important constraints on 
Skyrme-type effective interactions used to study neutron- rich systems. They do not, how- 
ever, provide information on the strength of the spin-orbit interaction, nor on other terms 
sensitive to density gradients, both of which may affect significantly the predicted properties 
of drip-line nuclei and of neutron-star matter. 

Ab initio calculations of finite nuclei, based on realistic models of nuclear forces, can pro- 
vide the necessary additional information, but they are more challenging. Recently [§J the 
energies of nuclei with A < 6 have been calculated essentially exactly with the Green's func- 
tion Monte Carlo (GFMC) method. Cluster variational Monte Carlo (CVMC) calculations 
have also been used to study 16 and the spin-orbit splitting (SOS) in 15 N Q. In this 
letter we report GFMC and CVMC calculations of states of seven and eight neutrons bound 
in a weak external potential well using the new Argonne two-nucleon |7j and Urbana three- 
nucleon interactions used in Ref. 0]. These interactions accurately reproduce the available 
two-nucleon scattering data and binding energies of A < 6 nuclei. Neutron matter is not 
bound, therefore an external well (V ex ) is necessary to hold the neutrons together. We have 
used a Woods-Saxon well with V Q = —20 MeV, R = 3 fm, and a = 0.65 fm, chosen such that 
with it alone only the single-neutron Is state is bound at -5.73 MeV, while the lp and higher 
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states are unbound. The investigated states of seven and eight neutrons are thus bound by 
both the well and the interaction between neutrons. We denote them by 8 n(J 7r = + ) and 
7 n(J n = 1/2- and 3/2-). 

The present variational Monte Carlo (VMC) and GFMC calculations are simpler than 
those for nuclei || because all nucleons are neutrons. The wave function is represented by a 
vector function of R (= fi,f 2 , ■ ■ ■ , r^) with 2 A spin components specifying the spin direction 
of each neutron. The VMC and GFMC calculations use a simpler variational wave function 
than those of Refs §,§]: 

*>, (1) 

Uij = u a (rij) Oi ■ <7j + u t {rij)Sij . (2) 

Here STI denotes a symmetrized product, Sy is the tensor operator, / c (r^) is the Jastrow 
correlation, and |<&) is an antisymmetric shell model wave function. The three-body cor- 
relations commonly used in nuclear tyy are omitted because they have little effect on the 
energies of low-density neutron systems, and the two-body spin-orbit correlations are dis- 
cussed later along with an improved \EV. The radial wave functions of the s- and p-orbitals 
in $ and the correlation functions f c , u a and u t are determined variationally. 

The GFMC calculations are carried out as described in Ref. Q with a simpler Hamilto- 
nian: 

H = - £ f-V* + E V„(i) + J2v' 8 (ij) + £ V ijk , (3) 

i i i<j i<j<k 

where the v' 8 does not contain L 2 or (L-S) 2 terms; it equals the charge-symmetric part of the 
Argonne fig interaction in the 1 5'o and 3 -Pj=o,i,2 two-neutron states. The small difference 
between the full vi$ and v 8 is treated as a first-order perturbation, whose contribution to 
the calculated energies is < 0.2 MeV. 

The calculated transient energies, E(t), 

E(t) = (^ v \He- iH - Eo)T \^ v )/(^ v \e~ {H - Eo)T \^v) , (4) 
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are shown in Fig.l. The E(r — > oo) converges to the lowest eigenvalue of the chosen J 71 ". In 
Fermi systems, the statistical error in E(t) increases with r as configurations diffuse across 
nodal surfaces of the wave function ||. Due to the large number of nodal surfaces in the 
wave functions of 7 n and 8 n states, it is difficult to study their transient energy for values of 
t > 0.04 MeV" 1 . The average values of E(t) for r = 0.032, 0.036 and 0.04 MeV" 1 , denoted 
by E, are shown by horizontal lines in Fig. 1. The E{r) of the 8 n(0 + ) and 7 n(l/2~) states do 
not have much r-dependence for r > 0.015 MeV -1 , suggesting that their E can be identified 
with the eigenvalues. In contrast, the E(t) of the 7 n(3/2~) state has more r-dependence. 
Consequently, the eigenvalue of the lowest 7 n(3/2~) state could be a little below its E value. 
However, we will neglect that difference and regard it as our best estimate of the eigenvalue. 
The GFMC estimate of the density distribution of neutrons in the 8 n(0 + ) drop is shown in 
Fig. 2. These results can be used to test the accuracy of the CVMC method and to further 
constrain the Skyrme type energy-density functionals used to study neutron-rich nuclei and 
neutron star crusts as discussed below. 

The CVMC method and its modification for SOS are described in Refs. The 
present CVMC calculations are more accurate; they include contributions of all correlations 
and interactions up to five-body clusters. In contrast, in |5|,|5| contributions of only static 
correlations and interactions were calculated up to four-body clusters and the momentum- 
dependent terms were evaluated only at the two-body level. With the simpler \Ev given by 
Eq. ([[]), the 1- to 5-body cluster contributions to the energy of 8 n(0 + ) state are respectively 
12.9, -54.5, 11.1, -3.8, and 1.1 MeV, which sum up to -33.3(2) MeV. The E(r = 0) is nothing 
but the variational energy calculated to all orders without cluster expansion. Its value of 
-33.7(1) MeV is very close to the CVMC result retaining up to five-body clusters. Note that 
even in this rather low-density system, the cluster expansion has a slow convergence and it 
appears necessary to include five-body cluster contributions to reduce the truncation error 
to < 2%. 

The simple (Eq. [TJ) is not very accurate; the energy obtained with it is ~ 4 MeV 
(or ~ 11%) too large. In CVMC we use the more general variational wave functions of the 
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where the three-body correlations, Uijk, are of the kind used in Refs. [|5]|| (note, however, 
that the commutator term is zero in pure neutron systems), and as before the Uij contains 
spin and tensor terms. The energies obtained with the ^>' v variational wave functions are 
respectively -35.6(1), -31.2(1), -29.7(1) MeV for the 8 n(0 + ), 7 n(±~) and 7 n(§~) states. They 
are only ~ 4 % above the GFMC energies -37.6(3), -32.3(2) and -31.2(2) MeV. The CVMC 
calculations require about a factor of 25 less computer time than the GFMC, even allow- 
ing for the variational search. It is not difficult to reduce the statistical error in CVMC 
calculations to a fraction of one percent. Much of the improvement in ty' v comes from the 
spin-orbit correlations omitted in the simpler \lV- I n the present GFMC calculations, the 
spin-orbit correlations are built in exactly via the propagation in imaginary time. 

Fragmentation of the P3/2 strength in 15 N results in the SOS of pi/2 and P3/2 quasi-hole 
states being ~0.6 MeV larger than the observed splitting between the lowest 3/2~ and l/2~ 
states in 15 N [§,|T0[. If the 3/2" hole strength in 7 n is similarly fragmented, the difference in 
the energies of lowest 7 n(l/2~) and 7 n(3/2~) states could be smaller than the SOS in 7 n. If 
W denotes the energy width of the fragmentation, the GFMC transient energies, E(r), for 
r > 1/W will include fragmentation effects. Since W is expected to be only a few MeV, it 
is very unlikely that the present E evaluated up to r ~ 0.04 has any fragmentation effects. 
The 7 n variational wave functions are constructed by removing an appropriate state from 
the |$) and thus correspond to quasi-hole states. Presumably, they too do not contain any 
fragmentation effects. Hence, we identify the difference between the calculated energies of 
the 7 n(l/2~) and (3/2~) states as the spin-orbit splitting. Its value is 1.1±0.3 and 1.4±0.1 
MeV in the GFMC and CVMC calculations respectively. 

The energies of uniform neutron matter calculated from several realistic models of nu- 
clear forces were plotted in Ref. [Q and compared with the results of four energy-density 



junctionals (EDF) used for astrophysical investigations involving dense matter: (i) Skyrme V 



- Vautherin-Brink Skyrme model 1 [O modified |I2] to fit the neutron- matter E(p) of Ref. 



[p"3|| ; (ii) SkM - Skyrme model M [rjj; (iii) FPS - a generalized Skyrme model fitted ap- 
proximately [|T5|,|I| to the nuclear- and neutron-matter energies of Ref. ||16|| ; and (iv) FPS21 
- a generalized Skyrme model [|I[] fitted accurately to results of Ref. [IB|. These EDF's re- 
produce the ground-state energies of stable closed-shell nuclei rather accurately. The root 
mean square deviations \AE/E\ between their prediction and experiment for 16 0, 40 Ca, 
48 Ca, 56 Ni, 90 Zr, 114 Sn, 140 Ce, and 208 Pb are listed in Table I. This table also contains their 
predictions for the energy of the 8 n ground state and the 7 n and 15 N SOS. The results for 
the density distributions p(r) of the 8 n(0 + ) state are compared with the GFMC and CVMC 
p(r) in Fig. 2. 

In trying to learn from the departures of the EDF results for 8 n and 7 n from our bench- 
marks, we concentrate on the FPS21 effective interaction, since it gives the closest fit to the 
neutron-matter energies. Possible sources of difference include the fact that the neutron- 
matter energies used (TB| date from an earlier period, whereas the benchmark results use 



newer interactions and more subtle computational techniques; also, and more importantly, 
the density-gradient terms in FPS21 are related to the effective-mass results [JTTJ assuming 
a zero-range nucleon-nucleon interaction with no spin exchange. The latter simplification 
is common to all of the EDF's we consider and is not well justified. On the assumption 
that this is a contributor to the discrepancy, we have examined the effect of adding a term, 
\a(p n y '(V 'pnf , to the FPS21 neutron EDF. The gradient term, \a{p n + p p ) /3 (Vp n - Vp p ) 2 , 
reduces to the form used for neutron drops, and would give little contribution for N ~ Z 
nuclei. However, in this work we have used it only for the neutron drops to avoid refitting 
the models to laboratory nuclei. Such a term can correct for the overbinding of 8 n(0 + ), and 
also reduce the central neutron densities. We find that it cannot give good fits to these quan- 
tities simultaneously, however. For the two sets of coefficients (/3, a) = (0, 150 MeV fm 8 ) 
and (2, 7 x 10 4 MeV fm 14 ) the ground state energies of 8 n are -39.6MeV and -40.2MeV 
respectively, and the neutron density distributions are shown as curves A and B in Fig. 2. 
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The present CVMC and GFMC results clearly indicate that the SOS predicted by the 
unadjusted Skyrme models for 7 n is too large, while it is good for 15 N. Relativistic mean-field 
models also predict a weaker spin-orbit potential in neutron-rich nuclei |[17|| . The neutron 
spin-orbit potential in these Skyrme models is of the form 

v£\r) = W*~(p(r)+p n {r)), (6) 

obtained by Vautherin and Brink JTlJ] assuming that it originates from two-nucleon inter- 
actions in the triplet-P state. Here p{r) = p n (r) + p p (r) and Wis is a constant determined 
from the SOS in laboratory nuclei like 15 N, and listed in Table I. The form of this potential, 
proportional to a radial derivative of the densities, indicates that apart from other depen- 
dences, the SOS obtained with a given EDF will depend on the nucleon central densities 
given by that model. The two modified versions of FPS21 just described, for the same value 
Wis = HOMeV used in Table I, each give a SOS in 7 n of 2.2 MeV. The reduction from 
the value of 3 MeV for unmodified FPS21 occurs because of the reduced central densities 
induced by the extra gradient term. The SOS in these modified models is still about double 
the value predicted by CVMC and GFMC, however. 

As discussed in Ref. 0, more than half the SOS in 15 N comes from three-nucleon contri- 
butions involving either two neutrons and a proton or vice- versa. In contrast, the three-body 
interaction and clusters give a very small contribution to the SOS in 7 n in CVMC calcula- 
tions. This suggests adding terms to the EDF that will produce a spin-orbit potential: 

VtHr) = W£ )l - Up (r) + p n {r)) + w[%- ±(p n (r)p p (r)) + W^} s - (p P (r)) 2 (7) 
r ar r ar r ar 

having separate two- and three-body contributions. The (p n (j)) 2 term is omitted because 
three-neutron clusters seem to give negligible V is ,(r). In neutron drops only the two-body 
Wff 1 contributes, while in N ~ Z nuclei, like 15 N, p n (r)p p (r) ~ p P (r) 2 and the sum W^ = 

(3) (3) 

Wi£ 8 + W^ls is the only relevant new parameter. 

With this modification to the spin-orbit interaction, and the FPS21 parameterization for 
the central interaction (including the gradient term for the neutron drops but not for 15 N), 
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we can fit the spin-orbit splittings of 15 N and 7 n exactly, using the parameter values W^ s ' 
= 61 MeV fm 5 and W^f = 745 MeV fm 8 . The spin-orbit splittings in the eight closed-shell 
nuclei mentioned earlier are modified only slightly a not unexpected result in view of their 
relatively small values of neutron excess 

In conclusion, we have made the first exact microscopic calculations of neutron drops 
in an external potential well. Our results suggest that the commonly used EDF's need 
modification in order to describe accurately neutron-rich nuclei. It appears that they predict 
neutron drops which are too dense and have too large a spin-orbit splitting. Additional 
density gradient terms need to be considered and the parameterization of the Skyrme 
must be modified to include three-body contributions. Our GFMC results also show that the 
CVMC using the improved \l/' y , with two- neutron spin-orbit correlations gives fairly accurate 
results. We plan to use CVMC to calculate the properties of larger neutron drops. This, 
together with data from stable nuclei, will provide a larger database for fitting a Skyrme 
EDF for studies of stable nuclei, neutron rich nuclei, and the surface of neutron stars. 
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TABLES 



TABLE I. Comparison of microscopic and Skyrme-model energies. 

8 n(0+) 7 n SOS Magic 15 N SOS a W? s 

Nuclei 

MeV MeV \AE/E\ rms % MeV MeV 5 



GFMC -37.6(3) 1.1(3) 

CVMC -35.5(1) 1.4(1) 6.1 b 

SkM -47.4 3.0 1.1 6.3 130 

FPS-21 -42.2 3.0 1.1 6.7 110 

Skyrme 1' -38.7 2.9 1.8 6.9 120 

FPS -32.5 3.5 1.2 6.7 110 



Experimental value 6.9 MeV deduced in Ref. ||. 
b With Argonne vi 4 , Ref. ||. 
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FIGURES 

FIG. 1. The transient energy E(t) for the GFMC calculations of 8 n(0+) and 7 n(l/2~ and 
3/2~) as a function of the imaginary time r. 

FIG. 2. Neutron density distribution for 8 n, according to methods described in the text. The 
curves labeled A and B come from modified versions of the effective interaction FPS21. 
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